Astronomy & Astrophysics manuscript no. NVSS dipole 


©ESO 2013 


January 28, 2013 





m 
o 

a 
in 

(N 



O 

u 

Oh 
i 

o 

o3 



(N 
> 
On 

in 
in 
in 

o 

m 

> 

X 

c3 



Cosmic radio dipole from NVSS and WENSS 

Matthias Rubart*, Dominik J. Schwarz** 

Fakultat fur Physik, Universitat Bielefeld, Postfach 100131, 33501 Bielefeld, Germany 
Preprint online version: January 28, 2013 

ABSTRACT 

We use linear estimators to determine the magnitude and direction of the cosmic radio dipole from the NRAO VLA Sky Survey 
(NVSS) and the Westerbork Northern Sky Survey (WENSS). We show that special attention has to be given to the issues of bias 
due to shot noise, incomplete sky coverage and masking of the Milky Way. We compare several different estimators and show that 
conflicting claims in the literature can be attributed to the use of different estimators. We find that the NVSS and WENSS estimates 
of the cosmic radio dipole are consistent with each other and with the direction of the cosmic microwave background (CMB) dipole. 
We find from the NVSS a dipole amplitude of (1.6 ± 0.6) x 10~ 2 in direction (RA,dec) = (154° ± 21°, -2° ± 21°). This amplitude 
exceeds the one expected from the CMB by a factor of about 3 and is inconsistent with the assumption of a pure kinetic origin of the 
radio dipole at 99.5% CL. 

Key words, observational cosmology, large scale structure, radio surveys, peculiar motion 



1. Introduction 

The assumed isotropy and homogeneity of the Universe at large 
scales is fundamental to modern cosmology. The isotropy is best 
seen in the cosmic microwave background (CMB) radiation and 
holds at the per cent level. The most prominent anisotropy of 
the CMB temperature is a dipole signal of AT/T as 10~ 3 . It is 
commonly assumed that this dipole is largely caused by the mo- 



tion of the Solar system through the Universe ( Stewart & Sciama 
1967). This interpretation seems to be fully consistent with the 



concordance model of cosmology. 

However, the observation of the microwave sky is not enough 
to tell the difference between a motion induced CMB dipole and 
dipole contributions form other physical phenomena, i.e. 



= d„ 



+ d 



primordial 



'ISW 



+ d 



foregrounds 



+ d n 



(1) 



Usually it is assumed that the primordial and the integrated 
Sachs-Wolfe (ISW) contribution to the CMB dipole are negligi- 
bly small and that foregrounds (the Milky Way) are under con- 
trol. Within the concordance model we expect a primordial con- 
tribution of ^primordial as 2 x 10~ 5 . The ISW contribution could 
be as large as 10~ 4 from the gravitational potentials induced by 
local 100 Mpc sized structures, without being in conflict with 
the concordance model (Rakic et al. 2006; Francis & Peacock] 
|2010) l. The noise term can be ignored due to excellent statistics 
of full sky observations. Thus the measured d cm ^ is directly used 
to infer the velocity o f the Solar system w .r.t. the CMB to be 
v = 369 + 0.9 km s _1 ( jHinshaw et al.|2009] >. It is used in many 
cosmological studies done in the CMB rest frame, e.g. supernova 
Hubble diagrams or measurements of large scale bulk flows. 

The effects of motion are not limited to the CMB, but should 
actually be detectable at any frequency. In order to test the hy- 
pothesis dc m b = rfmotiom it would be very interesting to measure 
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the dipole of another cosmic probe, such as that obtained by ra- 
dio point source catalogues. In this case one expects to find 



^radio — ^motion ^structure ^foregrounds d n 



(2) 



Besides the signal from our proper motion, we expect a signal 
from structure in the Universe and we expect a random dipole 
from Poisson noise. The dipole from structure is expected to 
dominate any catalogue limited to redshift z« 1. Thus we are 
interested in surveys with a mean redshift of order unity and a 
large enough sky coverage to be sensitive to the dipole. This 
makes radio catalogues the preferred probe to look at. Within 
the concordance model, the dipole signal induced by the large 
scale structure is then a subdominant contribution, as it is for 
the CMB. If we had a large enough catalogue, we could com- 
pare rf ra dio to d cm b. Any statistically significant deviation would 
be exciting, while finding a match would put the concordance 
model on firmer grounds. 

A first attempt to measure the radio dipole was performed 
by Baleisis et al. ( 1998 1 using a combination of the Green 
Bank ! 987 and the Pa rkes -MIT-NRAO catalogue s.|B lake & Wall| 
( [2T)Ml l, [Singal1 ( |2TJTT) and |Gibelyou & Huterer| ( |2012| i attempted 
to determine the dipole vector in the NRAO VLA Sky Survey 
(NVSS), with different conclusions. |Blake & Walt] ( |2002] l found 
a result that is in agreement with a purely kinetic origin of the 



cosmic radio dipole, but this was challenged by Singal (2011 



who finds a dipole amplitude four times larger than expected, 
but strangely enough pointing in a direction consistent with the 



CMB dipole. The analysis of |Gibelyou & Huterer| ( |2012[ ) finds 
both a different direction and an amplitude six times as large as 
the expected one. While Blake & Wall (2002i used a quadratic 
estimator, Singal ( 201 l| l and Gibelyou & Huterer ( |2012| l used 
different linear estimators to find the dipole direction. 

The purpose of this work is to discuss the use of linear es- 
timators of the cosmic radio dipole and apply several versions 
of them on the NVSS ( |Condon et al.|2002|l and the We sterbork 
Northern Sky Survey (WENSS) ( jRengehnk et al.|1997) . We re- 
solve the conflicts in the literature and extend the analysis to 
other linear estimators. 
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The NVSS survey covers about 10.3 sr of the sky and con- 
tains about 2 x 10 5 sources per steradian. For this survey the Very 
Large Array (VLA) in New Mexico (USA) has been used mea- 
suring at a frequency of 1 .4 GHz. The survey includes over 80 
per cent of the sky, missing only areas with declination 5 < -40°. 
The lower flux limit lies at 2.5 mJy for the 5<x detection of point 
sources. The NVSS was conducted by means of two different 
configurations of the VLA above and below 5 = -10°. 

The Westerbork Synthesis Radio Telescope in the 
Netherlands was operated at a frequency of 325 MHz to 
record the WENSS survey covering about 2.9 sr of the nothern 
sky and containing about 2.3 x 10 5 sources in total. This survey 
is made up of a main catalog for 6 e (28°, 76°) and a polar 
catalog for 6 > 72°. The 5cr detection limit for this survey is 18 
mJy. 

In order to analyze these surveys, we focus on linear esti- 
mators in this work. We do so for two reasons. Firstly, recent 
controversial results used linear estimators for the dipole direc- 
tion (Singal 2011; Gibelyou & Huterer 2012) and in one work 



also for the dipole amplitudes (|Singal|2011 1. Secondly, linear 



estimators are conceptually simpler. However, it is not expected 
that they are optimal (unbiased and minimal variance). The lin- 
ear estimators used in our analysis are asymptotically unbiased 
and their variance can be easily understood by analytic calcula- 
tions and by Monte Carlo simulations. 

The paper is organized as follows: First we discuss the ex- 
pected kinetic radio dipole. In section [3] we outline previous es- 
timates of the radio dipole. Linear estimators for full sky surveys 
are investigated in section [4j followed by a detailed analysis of 
the effects of incomplete sky coverage and masking in the next 
section. In section [6] we discuss the expected dipole amplitude 
from a flux based estimator. Our estimate of the radio dipole can 
be found in section|7]and is followed by a comparison with pre- 
vious results. We conclude in section[9] 



2. Kinetic Radio Dipole 

2.1. Doppler shift and aberration 



Ellis & Baldwin| ( [T98 4) predicted the kinetic contribution to the 
cosmic radio dipole for an isotropic and homogeneous cosmol- 
ogy. At redshift of order unity and beyond, we expect this kinetic 
contribution to be the dominant one. 

The spectrum of radio sources is assumed to be described by 
a power law, 



(3) 



where S denotes the flux and / the frequency. For simplicity, a 
mean value for the spectral index a is assumed. 

The number of observed radio sources per steradian depends 
on the lower flux limit and can be approximated by a power law 



dN . 
dH 



{>S)ocS~ 



(4) 



The value of x can be different for each survey. Typically x is 
assumed to be about one. 

Two effects have to be taken into account. The emitted radio 
frequency / rest is observed at the Doppler shifted frequency / D b S . 
The magnitude of this change depends on the angle 9 between 
the direction to the source and the direction of our motion, with 
velocity v. Observed and rest frame frequencies are related by 



where 5 is given by 
1 + v - cos(<9) 



5{v, 9) 



(6) 



Thus the observed flux changes due to our motion, since it de- 
pends on the frequency 



W/obs) oc 8f£ oc 6 1+a f£ <* S rest (/ obs )<5 



l+a 



(7) 



The first factor of 6 is due to the fact that the energy of an ob- 
served photon is enhanced due to the Doppler effect. 

Thus, the Doppler effect will change the number of observed 
sources above a given flux limit like 



ld_N\ _ld_N\ 
\d£l} obs U^Lt 



(8) 



Since the velocity of light is finite, aberration will also mod- 
ify the number counts. The position of each source is changed to- 
wards the direction of motion. The new angle ff (observed from 
Earth) between the position of the source and the direction of 
motion is given by 



sin# 



tan = 



■ cos 9 

Therefore, at first order in v/c, dQ. transforms like 

\2\ 



da' = dan 



2-cos6) + O 
c 



(9) 



(10) 



(11) 



fobs = /rest£(v, 9), 



(5) 



This can be combined with the Doppler effect to give the ob- 
served number density. After approximating S(v, 9) to first order 
in -, the result becomes 

c' 

i =(!) [i+p + *i + «)]PW> 

flfiobs \fl£2/ rest L \cJ 
The amplitude of the kinetic radio dipole is then given by 

d = [2 + x(l + a)](^\. (12) 

The kinetic radio dipole points towards the direction of our pe- 
culiar motion, which in an isotropic and homogeneous Universe 
must also agree with the direction defined by the CMB dipole. 

2.2. Expected kinetic radio dipole 

The measured CMB dipole is AT = 3.355 ± 0.008 tnK in the 
direction ( l,b) = (263.99° + 0.14°, 48.26° + 0.03°) dHinshaw| 
et al. 2009| l. In equatorial coordinates (epoch J2000) its direction 
reads (RA,dec) = (168°, -7° ). Compared to the CM B tempera- 
ture of T Q = 2.725+0.001 K ( |Fixsen & Mather|2002| l. this corre- 
sponds to a relative fluctuation of AF/T = (1.231+0.003)xl0 3 
and thus the velocity of the Solar system has bee n inferred from 
the C MB dipole to be v = 369.0 + 0.9 km s _1 ( |Hinshaw et~aL 
|2009l ). 

In order to find the expected amplitude of the kinetic radio 
dipole, we also need estimates for x and a. The typically as- 
sumed values are x = 1 and a = 0.75, which gives together with 
V = 370 km s~' a radio dipole amplitude of d = 0.46 x 10~ 2 . 
However, we can improve on that as x can be measured with help 
of the radio survey. Therefore we need to plot A^(> S ) against S 
like in figure [T] 
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Fig. 1. Number counts of the NVSS and WENSS surveys. A 
function f(S) oc S~ x is fitted to both data sets in the range of 
25 mjy < S < 200 mjy. Resulting values of x are 1.10 + 0.02 
for the NVSS survey and 0.80 + 0.02 for the WENSS survey. 



For the purpose of this work we find jcnvss = 1.10 ± 0.02 
and -*wenss = 0-80 ± 0.02. The mean spectral index cannot be 
inferred from the catalogues, as they provide data at a single 
frequency band only. We thus stick to a — 0.75, but include in 
the dipole error an uncertainty of Aa = 0.25 ( |Garn et al.|2 008). 
This results in the expectations: 



d NVSS = (°- 48 ± - 04 ) X 10 ~ 2 
^exp _ /n /it , nnn„ in-2 



(13) 

W1IN , S = (0.42 ± 0.03) x 10" z . (14) 
The error is dominated by the uncertainty in the spectral index. 



3. Previous results 

The first measurement of the radio dipole using the NVSS cat- 
alogue was performed by |Blake & Wall| ( |2002) >. In order to re- 
move corruption by local structure, all sources within 15° vicin- 
ity of the galactic disk have been removed. Additionally the 
clustering dipole contribution was reduced by ignoring sources 
within 30" of nearby known galaxies. The spherical harmonic 
coefficients a* s from the remaining NVSS catalogue have been 
determined up to I - 3. A model for a dipole distribution with 
an isotropic background has been constructed (aoo and ajo). Due 
to masking, this dipole distribution also influences higher multi- 
poles. After applying the same mask as for the NVSS catalogue, 
one finds a)™ del up to / = 3. A quadratic estimator (chi square) 
was used to compare the model with the observed coefficients. 

The resulting best-fit dipoles can be seen in table [T] The re- 
sults of Blake & Wall| p002) indicate a higher radio dipole than 
expected, however without statistical significance. 



Singal (201 1 1 used a linear estimator, originally proposed by 
|Crawford| ( |2009F , 



R 



3D 



and a variation of it, which we discuss below. For a large number 
of sources the isotropic background will clear away. The remain- 
ing vector /?3d will point towards the main anisotropy in the dis- 
tribution of number density over the sky. To get the correct dipole 
amplitude d one has to normalize this estimator depending on the 
number of sources. In Singal's analysis sources within 10° of the 



Flux 


N 


RA 


dec 


a 




(mJy) 




CO \ 

( ) 




o 


(10 -2 ) 


> 40 


125,603 


149 + 


49 


-45 ± 38 


0.7 ± 0.5 


1.02 




143 524 


161 ± 


44 


.77 + w 


+ 04 


74 


> 30 


166^694 


156 ± 


32 


2 ±33 


1.1 ± 0.4 


1.01 


> 25 


197,998 


158 ± 


30 


-4 ±34 


1.1 ± 0.3 


1.01 


> 20 


242,710 


153 ± 


27 


-3 ±29 


1.1 ± 0.3 


1.32 


> 15 


311,037 


148 ± 


29 


31 ±31 


0.8 ± 0.3 


1.81 


> 10 


431,990 


132 + 


29 


65 ± 19 


0.5 ± 0.2 


4.96 



Table 1. Best-fit dipole parameters from |B lake & W all (2002). 
Coordinate system and amplitude definition are adjusted for 
comparison with Singal's results (see tables|2]and|3]l. N denotes 
the number of point sources with flux above the indicated limit. 



galactic plane have been removed. In order to avoid a directional 
bias (see the more detailed discussion below) he reestablished a 
north-south symmetry of the NVSS by cutting all sources with 
dec > 40°. The results of Singal (201 1) are shown in table[2] The 
errors of the directional measurements are quite small here. This 
is an effect of an unexpectedly large amplitude, which simplifies 
the measurement. While the direction agrees with the one found 



by |Blake & Wall ( 2002| l, the dipole amplitude seems to be a fac- 
tor of about four higher than expected from the CMB dipole and 
twice as big as found by Blake & Wall (2002 1. 

Masking the supergalactic plane in order to reduce the con- 
tribution of local structure did not resolve the discrepancy. Since 
unknown clustering further away from the super galactic plane 
could also have contributed to the measurement, a second test 
was performed. A clustering contribution to the dipole would not 
give a signal proportional to cos 6. On the other hand, the differ- 
ence in number counts of areas that are opposite to each other 
should decrease with cos 9 (where 6 is the angle between an area 
and the measured dipole direction), if the measured dipole is due 
to our velocity. Singal was able to fit such a behaviour to the data. 
Therefore he concludes that the radio dipole amplitude is not due 
to local clustering . 
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N 


RA 




dec 
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o 
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171 ± 
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-18 ± 


14 


2.1 ±0.5 


> 40 


115,837 


158 ± 


12 


-19± 


12 


1.8 ±0.4 


> 35 


132,930 


157 ± 


11 


-12± 


11 


1.9 ±0.4 


> 30 


154, 996 


156 ± 


11 


-02 ± 


10 


2.0 ± 0.4 


> 25 


185,474 


158 ± 


10 


-02 + 


10 


1.8 ±0.4 


> 20 


229, 365 


153 ± 


10 


+02 ± 


10 


1.8 ±0.3 


> 15 


298,048 


149 ± 


09 


+15 ± 


09 


1.6 ±0.3 



Table 2. Dipole direction and amplitude from the number count 
estimator ( 15 1 from Singal (2011). 



Singal ( 201 \) also used a linear estimator for the distribu- 
tion of flux over the sky. This estimator is similar to the number 
density estimator ( 15 1, but weights each radio source by its flux 
Si, 



(15) /? flU x = £V<- 



(16) 



Like /?3 D , this estimator finds the main anisotropy and the am- 
plitude needs to be normalized. The brightest sources (S > 1000 
mJy) are removed, because they would dominate /?g ux other- 
wise. Results of this estimator are shown in table [3] The esti- 
mated directions are in agreement with the results of |Blake"& 
|Wall| ( |2002| l and the number count estimator results of Singal 
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( |201 However, the normalized dipole amplitudes d are even 
higher than those of the number count estimator R^. In section 
[6] we resolve this conflict. 



Flux 
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1000 > S > 50 
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11 


2.3 + 0.7 


1000 > S > 40 


1 14, 600 


159 ± 


12 


-11 ± 


11 


2.2 + 0.6 


1000 > 5 > 35 


131,691 


159 ± 


11 


-10± 


10 


2.2 + 0.6 


1000 > S > 30 


153,759 


159 ± 


11 


-07 ± 


10 


2.2 + 0.6 


1000 > S > 25 


184, 237 


159 ± 


10 


-07 ± 


09 


2.2 + 0.6 


1000 > S > 20 


228, 128 


158 ± 


10 


-06 + 


09 


2.1+0.5 


1000 > S > 15 


296,811 


157 ± 


09 


-03 + 


08 


2.0 + 0.5 



Table 3. Dipole direction and amplitude from the flux weighted 
number count estimator ( 16 1 from Singal (201 1 1. 



Most recently, Gibelyou & Huterer (201 2} measured a dipole 
amplitude (d = 2.7 + 0.5) X 10~ z towards (RA,dec) = (117 ± 
20°, 6 + 14°) from the NVSS. This direction is inconsistent with 
the studies mentioned above and the dipole amplitude is a fac- 
tor of five larger than expected. The authors used separate esti- 
mators for the direction and the amplitude. Their direction esti- 



mate is based on a linear estimator, originally proposed by Hirata 
p009] >, 



R 



3 DM 



ZNd V ~ 



(17) 



This three-dimensional estimator (3DM) is intended to be 
unbiased for arbitrary survey geometries and arbitrary masking. 
The idea is to achieve that with help of the second sum, which 
goes over Nr randomly distributed points, subject to the same 
masking. Therefore, the authors include all sources of the NVSS 
survey, except for those within 10° of the galactic plane. Below 
we show that this estimator has a direction bias, which depends 
on the real dipole anisotropy. 

We summarize, there is no agreement on the amplitude and 
direction of the cosmic radio dipole so far. 

4. Linear estimators for a full sky 

Let us first show that the estimator ( fT~5] > provides an unbiased 
estimate of the dipole direction. 

Starting from the distribution of the number of radio sources 
per solid angle (Hi, as seen by a moving observer in an oth- 



erwise isotropic Universe, the probability to find a given radio 
source within a solid angle dQ of position f is given by 

p(r)d£l= — (1 + rd)d£l, (18) 
4n 

where d denotes the dipole vector. 

In order to study the bias of an estimator we calculate its ex- 
pectation value with respect to an ensemble average. We do so 
below by means of Monte Carlo studies. For analytic consider- 
ations, for large N we replace the ensemble average by a spatial 
average, i.e. 



(1) 



r N 

= J [J dO^r,)! = 1, 



(19) 



thus we assume ergodicity. Note that the average is a linear op- 
erator. 



Now the expectation value of Crawford's estimator can be 
evaluated for large N, 



N N NT 

(Rid) = (J p <> = X< ? <> = ^ J dfl (1 + r ■ d) i 



(20) 



This calculation holds for independent, identically distributed 
positions r,, thus without clustering effects. Only the second term 
survives the integration and thus the expected dipole estimator is 



(R 3D ) = ^Nd. 



(21) 



Naively, one could now estimate the dipole signal by d^u = 



|*3D. 



We conclude that d^o provides us with an unbiased estimate 
of the dipole direction d for a full sky sample. However the esti- 
mated dipole amplitude \dso\ is biased. 

To understand the origin of this bias let us first consider 



<dj D 



>= [1- -)d 2 + - >d 2 . 
N N 



(22) 

The inequality holds for large N and d < 3 (in case of large 
dipole amplitudes [d = (9(1)] our ansatz (19i should also take 
many-point correlations into account). Thus d^ D is definitely bi- 
ased towards higher amplitudes. However, in order to prove that 
\dio\ is biased, we would need to calculate (I^dD- We do this by 
means of the random walk/flight method. 



4. 1 . Random flight 

Adding up vectors for each point of a survey corresponds to a 
random walk with unit step size. To be more precise this is a 
random flight, since the problem is three dimensional. Even for a 
vanishing dipole, such a random flight is unlikely to return to the 
origin after N steps. This describes the noise of any realisation 
of an isotropic distribution of N sources. 

Following Crawford ( 2009 1, we determine the distance r 
from the origin after N steps from the probability density of a 
random flight process 



P N (r)dr 



54 



nN 3 



1/2 



r 2 exp 



3r^ 
~2N 



dr. 



(23) 



The probability of measuring a dipole signal of an amplitude 
bigger than R in a random flight is 



Pn(R > Rpd) 



f 



drP N (r) = 1 - pCL. 



(24) 



A confidence level pCL can be choosen, leading to errorbars 
for a measured dipole vector R^ D + R p cl- To estimate the di- 
rectional uncertainties of this method, Crawford (2009) made 
the following argument: at a given confidence level the random 
flight corresponds to a step of length up to R p cl- Adding R p cl 
perpendicular to the measured dipole R^u allows us to estimate 
the maximal offset in direction. Using trigonometry, one can re- 
late Rpci to the directional uncertainties expressed as the angle 



80, 



R 



pCL 



arcsin ■ 



pCL 



Rm 



(25) 



The expected magnitude of the random dipole contribution 
is estimated from d23l as 



/ orandom \ 
<^3D > 



f 

JO 



54 



nN 3 



1/2 



r 2 exp 



3r* 
'2N 



dr * 0.92 V/V. (26) 
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Fig. 2. Ampliude bias of the full sky estimator dy,u- Data repre- 
sent mean and empirical variance of 1000 simulations for each 
N. A function d ohs (N) = yjD 2 + 9A 2 /N is fitted to the simulated 
data, with best-fit values A = 0.908+0.002, D = (0.451 ±0.001 )x 
10~ 2 . The expected dipole amplitude (d = 0.0046) is indicated 
by the horizontal line. 



telescopes are limited to observe at certain declination ranges. 
Additionally, the Milky Way is masking parts of the sky. Here 
we will discuss some of the effects caused by incomplete sky 
coverage. 

The upcoming Low Frequency Array (LOFAR) Tier-1 sur- 
vey will cover about half of the sky (2n), thus we will first focus 
on this situation. As a second step we generalize this to an ar- 
bitrary axisymmetric survey geometry and include the effects of 
masking. 

5. 1 . Random walk 

Let us assume a survey geometry that covers all of the Northern 
hemisphere and ask how the estimator of the radio dipole ( [T5| 
has to be modified. For the first two Cartesian components of 
/?3D there should be no systematic problem, but the third com- 
ponent will definitely be biased. It is necessary to remove the 
effect of the incomplete sky from this z component. Consider 
the expectation value of /?3d for the Northern hemisphere 



N N „ 

<^3D>hemi S phere = ( / = — — dCl(l + f ■ d)f, 

j~t 2 7T J<5>0 



(28) 



Since the random dipole has no distinguished direction, there 
is no direction bias of the linear estimator for a full sky map. 

Even for vanishing d, this gives rise to a non-vanishing esti- 
mate of the dipole amplitude t/™ ldom = 2.8N~ 1/2 and is thus the 
origin of the amplitude bias. 

Motivated by ( 22 1 we make the following ansatz for the 



dipole amplitude and its bias: 



<d 3D ) = \\d 2 + T ^<*SS dom > 2 . 



9 



(27) 



In order to verify this analytic estimate of the biased amplitude, 
we simulated full sky maps including a velocity dipole (with v = 
370 kins -1 ). From these simulated catalogues we extracted the 
observed amplitude d \, s . Figure [2] shows the simulated data, the 
true value of the dipole amplitude and a fit of the form f(N) = 

V£> 2 + 9A 2 /N. 

The best-fit values are A = 0.908 + 0.002 and D = (0.451 ± 
0.001) x 10~ 2 (statistical errors only). These numbers should be 



compared to the factor 0.92 from ( 26 1 and the simulation input 
of d = 0.46 X 10~ 2 . The reduced chi- square of the fit is 7 x 10~ 5 . 
Thus the Monte Carlo simulations agree with the theoretically 
motivated ansatz ( 27 1 for the expected dipole amplitude of the 
estimator d^u- 

We conclude that even for a perfect full sky catalogue (no 
masking, complete in flux, perfect flux and position measure- 
ments), the amplitude of the linear estimator is biased towards 
higher values. Increasing the number of radio sources will re- 
duce the bias. The estimator d^u is asymptotically unbiased, but 
this is of limited practical use for the analysis of NVSS and 
WENSS data. A similar bias of the dipole amplitude is also 
found for the other linear estimators introduced above. This find- 
ing is in agreement with Gib elyou & Huterer| ( |2012| l, who use a 
linear estimator to find the direction of the dipole. 



5. Linear estimators for an incomplete sky 

So far we assumed full coverage of the radio sky. More realis- 
tic catalogues cover just a fraction of the sky, as all earth based 



where p c [ = 1/[1 + (d/2) cos H-j)] accounts for the proper nor- 
malisation of the probability distribution on the hemisphere in 
presence of a dipole. An in ( p~8^ > becomes 2n for obvious reasons. 
The integral can be evaluated by hand. One finds 



(J?3D hemisphere - Npd 



' \d cos ipd sin ftj 
Ic/sin % sini?^ 



+ TffC0S?°V 



(29) 



where §d and tp& denote the dipole position in spherical coordi- 
nates. The z direction is strongly influenced by the incomplete 
sky coverage. The total number of observed sources is not inde- 
pendent from the amplitude and orientation of the dipole itself. 

Nevertheless there is no problem in the evaluation of <p^ be- 
cause one can calculate 



tpd = arc tan 



Ry 



(30) 



Here N as well as sin cancel out. So the 2D estimator is un- 
biased with respect to tp c /. Therefore we propose a pure two di- 
mensional estimator: 



R 



2D 



z 



cos tpi sin i?; \ 
sin <pi sin § t 




(31) 



From this one can still use (R 2 + R 2 y 1 ^ 2 for evaluating 
ds,m-& l jNp t i/3 and Let us take a look at the factor sini?^. 
Sources near the pole will make a smaller contribution than those 
further away. If a source near the pole is shifted by a small dis- 
tance, the value (/?, of this source could change dramatically. So 
the weighting terms compensate for this artificially big errors, 
which are just a relic of the coordinate system. 

Let us now estimate the uncertainties of the estimator R2d- 
The problem corresponds to an isotropic random walk process 
with variable step size. The probability density for a displace- 
ment of r for such a random walk is 



PN(r)iDdr = — rexp 



'2N 



dr. 



(32) 
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Similar to the random flight we determine R p cl defined by 

P N (R 2D > R p cl) = I drP N ( r ) 2D = 1 - P CL. (33) 

Here again pCL is the confidence level. It is possible to solve the 
above integral analytically 



f 



dr— r exp 

N 







f 3R Icl) 




= exp 


{ 2N J 



So R p cl is given by 



R 



2N 



pCL - 



In 



1 



1 -pCL 



(34) 



(35) 



Now one can use the same argument as for the random flight to 
evaluate the uncertainty of the ipd estimation 



fypCL = arcsin 



R 



pCL 



Rid 



(36) 



In this way one can directly determine error bars for measured 
results of ipd calculated via ( 30 1. Using this two dimensional esti- 
mator one cannot measure the dipole amplitude d itself but only 
the combination ds,m& d Npd/3. Therefore it can only give an 
lower limit for d. 

Like in the case of the full three dimensional estimator, this 
version is also biased in the measurement of the amplitude. The 
expectation of the random contribution can be calculated via 

:D ,,- , r-rexp|-— ]dr*0.72Vw. (37) 



/ orandomx 
\^2D > 



2N 



So we expect our estimator to measure a combination of this 
random contribution and the true velocity dipole and make the 
ansatz 



(J 2D ) = y \d 2 sin 2 & d + — (R™ dom ) 2 , 



(38) 



where we used pd « 1. Like above, we verify this via Monte 
Carlo simulations, shown in figure [3] 

5.2. Direction bias 

For any masked or incomplete map of the sky, we cannot 
measure the mean source density N/(4n), i.e. the monopole. 
Therefore we always have to keep in mind that the observed 
mean density is just an approximation. Gibelyou & Huterer 
( 2012| l have used an estimator proposed by Hirata ( 2009 ), which 
implicitly assumes, that the monopole is known. Based on the 
knowledge of the monopole, this estimator would compensate 
for masking effects and incomplete sky coverage by substract- 
ing a pure random isotropic map from the observed dipole term 
via ( 17 1. However, we cannot know N. 



For the previous example of a hemisphere, the monopole 
density is N/(4n) = Np d /(2n), where N = N D in (fl7}. The ran- 



dom map can be made up of an arbitrary number of sources Nr, 
but is reweighted by the observed number of sources Nd/Nr, 
instead of N/Nr. This introduces a directional bias, in addition 
to the previously discussed biased amplitude. One can see this 
explicitly by applying ( [17) ) on our model of a dipole modified 
isotropic hemisphere. We find 



ere 



N Pd 



cos ip d sin § d 
sin ifd sin § d 
cos $ d 



) 1 


( 1 


+ 2 N 







, Pd ~ 1 / 



(39) 
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Fig. 3. Ampliude bias for the estimator d 2 o on a hemisphere. 
Data represent mean and empirical variance of 1000 simulations 
for each N. A function d ohs (N) = tJd 2 + 9A 2 /N is fitted to the 
simulated data, with best fit values A = 0.712 + 0.003 and D = 
(0.444±0.002)xl0~ 2 . The dipole amplitude (0.0046) is indicated 
by the horizontal line, the dipole vector is assumed to lie in the 
equatorial plane (sini?^ = 1). 



which clearly is not parallel to d. For small values of d we can 
expand p d and obtain 



_ N 

(^3DM)hemisphere — l^d 



' cos tpd sin §d N 
sin ip d sin §d 
2 cos & d 



+ 0(d 2 ). 



(40) 



The z component of the dipole is underestimated by a factor of 
4 for the geometry of a hemisphere. Despite cancelation of the 
leading term of the bias of the z direction, the dipole direction 
remains biased. Less symmetric survey geometries lead to a bias 
of all dipole components for this estimator. 

The best strategy to avoid any directional bias for a linear es- 
timator is to make the survey geometry point symmetric around 
the observer for three dimensional estimators like /?3 D or point 
symmetric around the zenith in case of two dimensional estima- 
tors /?2d- This implies for the NVSS that one has to cut sym- 
metric in declination, such that both polar caps are missing, a 
strategy that was used by Singal (2 01 1\ . 



5.3. Masking 

The use of a masked sky additionally affects the dipole measure- 
ment. In general, the estimated dipole direction and amplitude 
depend on the position of the true dipole relative to the mask. 
Cutting areas with large dipole contribution will reduce the am- 
plitude and vice versa. 

In the following we consider masks that are point symmetric 
with respect to the observer for all considered 3-dimensional es- 
timators, respectively point symmetric with respect to the zenith 
for the 2-dimensional estimators. Constructed in such a way, 
masking does not introduce any directional bias. Nevertheless, 
we have to consider the effects of masking on the estimated 
dipole amplitude. 

A simple method to correct for this effect was put forward 
by Singal ( |201 l| l, who introduced a masking factor k, 



R 



3D 



(41) 
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Such a factor is expected to be a function of shape and position of 
the mask as well as the dipole position. In most cases an analytic 
calculation of this factor is impossible. For simple mask geome- 
tries some analytic results can be found in Rubart|( |2012| l. An al- 
ternative approach is to simulate the effects by means of Monte 
Carlos and to compare simulations with and without masking. 
The ratio of both results after a large number of simulations pro- 
vides the masking factor k^u- 

Doing so, we found some interesting effects. The masking 
factors depend on the number of objects in the simulated maps 
as well as on the true dipole amplitude. This can be explained as 
follows. The masking will mainly affect the kinetic dipole contri- 
bution, while the random dipole will increase due to the decrease 
of the number of objects in the masked catalogue. Therefore we 
expect the amplitude to be 



(^3D)mask - \l ^d 2 ' 



9 

Jf2 



</?™ dom > 2 . 



(42) 



£30 should only depend on the properties of the mask and not 
on the amplitude of the dipole. To test this, we created simulated 
maps for two different dipole amplitudes, with values motivated 
by the CMB measurement drum = 0.46 x 10~ 2 and by the mea- 
surement of |B lake & Wall| ( |2002| at 25 mJy d B w = 1.1 X 10~ 2 . 
We used the direction RA = 158°, dec = -4° in both cases. The 
mask agrees with the one used by Singal ( 201 1), i.e. we removed 
all sources with \6\ > 40° and \b\ < 10°. Resulting amplitudes for 
different numbers of sources are shown in figure |4] 
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Fig. 4. Amplitude bias of the 3-dimensional estimator for the 
masked NVSS geometry of Singal] ( [201 l| l. Data represent mean 
and empirical variance of 1000 simulations for each N. A func- 
tion d ohs (N) = yJ(KD) 2 + 9A 2 /N is fitted to the simulated data, 
with best-fit values A = 0.883 ± 0.006, KD = (0.642 + 0.005) x 
10~ 2 and A = 0.847 ± 0.016, KD = (1.59 ± 0.008) x 10~ 2 for 
the expected kinetic radio dipole and the radio dipole measured 
by Bl ake & Wall| ( |200"2| , respectively. The simulated dipole am- 
plitudes, without masking, are indicated by the horizontal lines. 



A similar argument holds for the two dimensional estimator. 
Here we expect a behaviour of the form: 



<</2D>mask = ^Afpd 2 SUV 2 & d + — <fl™ d ° m > 2 . 



(43) 



Using the same assumptions about the dipole term and the same 
mask as before, we obtain figure [5] 
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Fig. 5. Amplitude bias of the 2-dimensional estimator for the 
masked NVSS geometry of Singal (201 1 ). Data represent mean 
and empirical variance of 1000 simulations for each N. A func- 
tion d obs (N) = yJ(KD) 2 + 9A 2 /N is fitted to the simulated data, 
with best-fit values A = 0.810 ± 0.009, KD = (0.589 ± 0.008) x 
10~ 2 and A = 0.745+0.014, .KD = ( 1.493+ 0.006) x 10' 2 , for the 
expected kinetic radio dipole and the dipole measured by Blake 
|& Wall| ( |2002} , respectively. The simulated dipole amplitudes, 
without masking, are indicated by the horizontal lines. 



This time, we find &2d = 1.3 for both assumed velocities. 
The simulations support our assumption that the masking factor 
does not depend on the dipole magnitude d. 

5.3.1. Masking correction 

Although the masking factor k does not depend on the amplitude 
d, it may depend on the dipole direction d. Therefore it would be 
necessary to repeat the analysis of the previous section for each 
dipole direction found. In order to reduce the simulation effort, 
we rely on the following method instead. 

For the full, as well as for the masked sky, surveys with 10 6 
sources were simulated. This choice guarantees that we inves- 
tigate masking effects and not effects due to shot noise. The 
mean dipole amplitudes are determined for 10 3 simulated full 
and masked sky surveys, respectively. The ratio of the masked 
sky mean amplitude to the full sky mean amplitude is denoted k. 
This ratio provides a first approximation to the masking factor. 



First of all we can conclude that (42 1 is a good fit for the 



behaviour of the simulated maps in both cases. The measured 
amplitudes are larger than those estimated from full sky maps. 
We can now calculate £30- In both cases it turns out to be 1.4, 
which could be used to correct the amplitude estimate. 



k = 



^Jk 2 d 2 + ^</? random ) 2 
^jd 2 + j^{R™ niom ) 2 



(44) 



The influence of the random dipole tends to bias k towards 1 (as 
can be easily seen in the limit of a small number of sources). 
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This bias can be compensated by rewriting the above formula 
into 



k = Jk 2 + 



0,^random^2 



d 2 N 2 



(k 2 - 1). 



(45) 



The values of N and d are input parameters of the simulation. 
From the last section we know the values of <fl random > for the 
three- as well as for the two-dimensional case. Therefore we can 
transform the approximated masking factor k into the unbiased 
masking factor k. 

6. Dipole from a flux weighted estimator 

Let us now turn to a discussion of the flux weighted estimator 



16 1, which was used by Singal (20TTJ). 

We first take a closer look at the theoretically expected value 
of the flux dipole df\ ux . For simplicity we assume full sky cover- 
age. For a large number of sources, 



dS — —S cos 6, 
dQdS 



(46) 



where 9 is the angle between a source and the dipole direction d. 

We now determine the number of sources per flux and solid 
angle as a function of the observer velocity. At zeroth order in 
velocity, this density n(S) is isotropic. As for the number counts, 
stellar aberration and the Doppler effect have to be taken into 
account. Stellar aberration gives rise to 



d 2 N 
dOdS 



n(S)(l +2-cos#). 

c 



(47) 



The relativistic Doppler effect alters the observed fluxes. When 
we observe a source in the direction of motion, we measure a 
higher flux than if we were at rest with respect to the isotropic 
and homogeneous Universe. We relate the observed flux So to 
the flux that is measured by an observer with vanishing peculiar 
motion, 



S rest ~ S - S (l + Or)- COS I 

c 



(48) 



We assume the power law n(S) = aS x to hold for observers at 
rest and Taylor expand around the observed flux 



n(S KSt ) ~ n(S )[l + x(l + a)- cos 6]. 

c 



(49) 



If the assumed power law holds for all sources of a survey, then 
x = 1 + x, with x as defined above in Q. Combining the Doppler 
effect and the stellar aberration leads to 



d 2 N v (ivy 

rio(S )[l +(2 + x(l + a))- cos 0] + OU-j 



dCldS 



(50) 



This result only holds under the assumption of a power law be- 
haviour of the number counts. It is crucial to keep this in mind. 

6.1. NVSS 

Let us now see, if we are allowed to make this assumption for 
the analysis of NVSS data. The plot in figure [6] demonstrates 
that the power law is not valid for the flux range used in the 
analysis of Singal (2011), i.e. fluxes up to 1 Jy. Actually the slope 
steepens for the larger fluxes considered. The best fit power-law 



gives a reduced chi-squared value of x 2 = 122. That means that 
the observed data cannot be fitted by a power law. 

Thus Singal's assumption x « 1 does not hold for two rea- 
sons. Firstly, for a pure power law we would expect x = l+i * 2, 
which is close to what we find: x. m 1.9. Secondly, the power-law 
assumption only applies to about one half of the data. 
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Fig. 6. Differential number counts of the NVSS catalogue, 
S min - 10 mJy, best fit values for f(s) = a ■ s~ x with 25 mJy < 
S < 1000 mJy are a = 2.1 x 10 7 and x. = 1.9. 



In conclusion, the result of unexpectedly large amplitudes in 
Singal ( 201 1 ) could partly be explained by this two effects. With 
power law behaviour one should use x = 1.9 or larger to also 
account for the steepening of the spectral index at high fluxes, 
which increases the expectation value by 



<^obs)l.t= 



1.9 



<^obs)Lv= 



> 1.4. 



(51) 



The results of the flux weighted estimator from Singal ( 201 1\ 
should be reduced by at least a factor of 1 .4. 



Compared to the estimator of Crawford, the estimator ( 16 1 
stresses sources with high flux. To avoid the domination of 
a small number of sources, sources with S > 1000 mJy are 
not taken into account. Nevertheless, this estimator is stressing 
bright sources. These sources are, on average, nearer than the 
rest. Hence ( 16 1 might be dominated by nearby sources and by 
atypically bright ones. This seems to be yet another weakness of 
this estimator, since the local universe is anisotropic. 



7. Dipole estimates from NVSS and WENSS 

7.1. 3D linear estimates 

We are now in a position to check the three dimensional esti- 
mations of the radio dipole in the NVSS catalogue. As we have 
shown above, the estimator used by |Gibelyou & Huterer ( 2012| l 
( 17 1 gives rise to a biased dipole direction and thus is not fur- 
ther considered in this work. The flux weighted estimator ( |T6] > is 
also of limited use, as the NVSS data cannot be described by a 
power-law over all fluxes of interest. We thus focus here on the 
simplest linear estimator ( [15) , 

In order to obtain an unbiased direction estimate, the cut sky 
is adopted. The masking factor k is 



geometry of Singal (2011 



determined for every measured dipole anisotropy direction, as 
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described in |5.3.1 The pure estimator results c/30 are then cor- Counter Galactic plane (CG). A HEALPi^map of the remain- 
rected for the masking bias and we obtain d^g. ing NVSS sources can be seen in figure[7] The color of the pixels 

encodes the number of sources per pixel. 



Flux 
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RA 


dec 


d iD 
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(mJy) 
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(10- 2 ) 


50 


91,662 


170 ± 25 


-17 ±25 


2.53 


1.38 


1.8 ±0.8 


40 


115,917 


156 ±28 


-18 ±28 


2.03 


1.29 


1.6 ±0.8 


35 


133,022 


156 ±25 


-11 ±25 


2.24 


1.32 


1.7 ±0.7 


30 


155, 120 


156 ± 21 


-2 ±21 


2.39 


1.35 


1.8 ±0.7 


25 


185,649 


158 ± 21 


-2 ±21 


2.16 


1.34 


1.6 ±0.6 


20 


229, 557 


153 ± 19 


2± 19 


2.10 


1.30 


1.6 ±0.6 


15 


298,289 


149 ± 19 


17 ± 19 


1.84 


1.28 


1.4 ±0.5 



Mollweide view 



Table 4. Dipole direction and amplitude from NVSS. The esti- 
mator (JT3J was used. Excluded are sources with \S\ > 40° and 
\b\ < 10° (J2000). 



All results for right ascension agree within their error bars. 
The same holds true for the dipole amplitudes. This self- 
consistency indicates the absence of significant systematic er- 
rors. Unfortunately, we can not state the same for the declination 
results. One observes an significant increase in declination with 
respect to deceasing flux limits. This effect is very likely due to 
a relic of the NVSS survey procedure. Sources below 6 = -10° 
were measured by means of a different alignment of the Very 
Large Array. The source density is therefore smaller in this area. 
Therefore it is expected that the dipole measurement will show 
increasing values of declination at the fainter flux limits. This 
makes it hard to trust the declination results below a flux limit of 
25 mJy, (note that there is no significant difference between 20 
mJy and 30 mJy). 

The results of table [4] can be compared to those in table [2] 
Number of sources and direction results are almost the same. 
Deviations could be explained by minor differences in the pre- 
cise form of the mask. We used a different method for estimating 
uncertainties in the direction and amplitude measurement, than 
those used by Singal (201 1 ). Since our method (described above) 
is more conservative, we obtain larger errorbars. All dipole am- 
plitudes in this table are slightly below those from table [2] In 
Singal] ( |201 l| l a different method was used to obtain the dipole 
amplitude from the estimator (15) , which can explain this dis- 
crepancy. In principle we can recover the results of Singal ( 201 1 ) 
and confirm an unexpectedly high dipole amplitude. 



7.2. 2D linear estimates 

A disadvantage of the NVSS catalogue is that the sampling depth 
changes at S — -10° (less sensitivity at lower declinations). This 
could lead to a directional bias of the NVSS data analysis and 
thus it is interesting to also use the two dimensional estimator 
presented in this work, as this effect cannot lead to a bias in 
this case. For the WENSS analysis, a three dimensional linear 
estimator cannot avoid directional bias. 



7.2.1. NVSS 

A major advantage of the estimator Rod is that it does not re- 
quire a north-south symmetry of the catalogue. Therefore, the 
declination limit of the NVSS catalogue is no problem. 

As the estimator R 2d requires a point symmetry around the 
north pole, we cannot remove the Galactic plane only. For each 
removed point we also need to subtract the point which is 180° 
away. When we do so, a second plane occurs, which we call the 




Fig. 7. A map of the number counts in HEALPix pixels from 
NVSS. The pixel size corresponds to /V s ;d e = 32. Shown are 
equatorial coordinates at epoch J2000. The NVSS contains data 
at 6 > -40° and the galactic plane and a "counter galaxy" are 
masked (CG mask) in order to avoid galactic contamination and 
to restore point symmetry with respect to the zenith. 



An alternative would be the mask used by Singal] (2011 1. 
That mask is a combination of two great cycles and would there- 
fore also work for /?2d- However, it turns out that the mask with 
the CG removes fewer sources. Thus we decided to use the CG 
mask. 

The next step is to evaluate the masking correction k of the 
CG mask. As this factor depends on the right ascension and on 
the declination of the dipole anisotropy, we need some additional 
information. From /?2D we estimate the right ascension. As the 
declination cannot be evaluated with Rid, we use the declina- 
tion as provided by the three dimensional estimator in order to 
determine k. These values are certainly not exact, since a differ- 
ent mask is used now. The influence of a small change in dipole 
declination on the evaluated factor k is discussed in section lT.2.21 

We reduce the dipole amplitude of the 2D estimator by the 
masking factor k and obtain the masking corrected amplitude 
t/jD - The results of this procedure for the NVSS catalogue and 
different flux limits can be found in table [5] 

Obtained right ascensions and amplitudes are stable with re- 
spect to different flux limits. For all flux limits > 15 mJy, the esti- 
mated right ascension of the radio dipole is in agreement with the 
CMB prediction of RA = 168°. Only when we include sources 
as faint as 10 mJy, find a 3cr deviation. However, we know that 
the catalogue is incomplete at its faint end. 

The masking corrected dipole amplitudes d^ are signifi- 
cantly above the CMB prediction. To some extent this is ex- 
pected due to the discussed amplitude bias. A detailed discussion 
is presented in the next section. 

7.2.2. WENSS 

We finally present the first estimation of the radio dipole from the 
WENSS catalogue. We cannot use the three dimensional linear 



http://healpix.jpl.nasa.gov 
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Table 5. Dipole right ascension and amplitude dsindj from 
NVSS. The 2D estimator pTj ) is used and all sources with 
5 > -40°, except the galactic and counter galactic planes (CG 
mask), are included. 



estimators on this catalogue, since it only contains sources with 
dec > 28°. The two dimensional estimator on the other hand can 
be used here. Again, we need to remove the Galactic plane and 
the Counter Galactic plane to reestablish a symmetry around the 
north pole. The remaining WENSS catalogue is shown in figure 



8. 
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Fig. 9. Differential number counts of the WENSS catalogue, 
Smin - 5 mJy, best fit values for f(s) - a ■ s~ x with 25 mJy < 
S < 1000 mJy are a = 2.6 x 10 6 and x = 1.6. 



is relatively small, compared to shot noise uncertainties. We as- 
sume dec = 0° for the determination of k, based on the expecta- 
tion from the CMB dipole and the NVSS radio dipole estimates. 
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Table 6. Masking correction k for WENSS with CG mask and a 
dipole with/M = 120°. 



Fig. 8. A map of the number counts in HEALPix pixels from 
WENSS. The pixel size corresponds to N s a e = 32. Shown are 
equatorial coordinates at epoch B1950. The WENSS contains 
data at 6 > 30° and the galactic plane and a "counter galaxy" are 
masked (CG mask) in order to avoid galactic contamination and 
to restore point symmetry with respect to the zenith. 



In order to choose the best flux limit, we analyse the differ- 
ential number counts between 10 mJy and 1000 mJy. Figure [9] 
shows that the WENSS catalogue seems to be incomplete be- 
low 25 mJy. In Rengelink et al. (1997]) the completeness of the 
WENSS catalogue is claimed to hold only above a limit of 30 
mJy. From figure|9]we infer that a simple power law also applies 
to the source counts between 25 and 30 mJy and thus we include 
sources down to a flux limit of 25 mJy in our analysis. 

We cannot obtain information on the declination of the 
dipole from the WENSS catalogue by means of the two di- 
mensional estimator applied in this work. This could in prin- 
ciple be a problem for the determination of the masking factor k. 
Therefore, we further investigated the effect of different dipole 
declinations. Assuming the WENSS symmetry and a right as- 
cension of 120° (close to the results given in table [7]), we calcu- 
lated k for 7 different values of declination, see Table [6] For this 
mask, the dependence of k on the right ascension of the dipole 
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(io- 2 ) 


40 


67,052 


124 + 51 


1.31 


0.45 


2.9 + 2.3 


35 


73, 653 


123 + 46 


1.36 


0.47 


2.9 + 2.1 


30 


81,863 


122 ± 48 


1.24 


0.45 


2.8 + 2.1 


25 


92, 600 


117 + 40 


1.36 


0.47 


2.9+1.9 



Table 7. Dipole estimate from WENSS based on 2D estima- 
tor using peak flux values for all sources with 6 > 30°, except 
those in the galactic and counter galactic planes (CG mask). Our 
WENSS analysis uses positions at epoch B1950. 



The results of the WENSS analysis are presented in Table 7] 
Although the WENSS catalogue covers only one fourth of the 
sky, we find that it can be used for the estimation of the radio 
dipole. A problem is the limited number of sources that are left 
after masking the galaxy and restoring the required symmetry 
of the catalogue. This leads to larger error bars, compared with 
the NVSS estimates. We can conclude that the observed dipole 
anisotropy in the WENSS catalogue is in agreement with the 
measurements from NVSS, which is a nontrivial statement, as 
we are probing radio sources at different frequencies (325 MHz 
vs. 1.4 GHz). 
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simulation i i 



0.005 0.01 0.015 0.02 0.025 
d 



0.09 
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0.01 




simulation i i 



0.005 0.01 0.015 0.02 0.025 
d 



Fig. 10. Histogram of dipole amplitudes for 100,000 simulations of the three dimensional (left) and two dimensional (right) esti- 
mator, assuming the CMB expectation and a slope of x = 1.1, with 185,649 (left) and 195,245 (right) sources per simulation and 
appropriate masks. The black vertical lines are the NVSS results of this work. 



8. Comparison of results 

We summarize the various results from the literature and this 
work in table [8] The results of this work are highlighted with 
bold faced letters. For comparison we focus on the flux limits of 
25 mJy and 15 mJy. 

All estimated dipole directions, both from the NVSS and 
from WENSS are in good agreement with each other and with 
the direction from the CMB dipole, with the exception of the re- 
sult from Gibelyou & Huterer ( 2012| l. As explained in section 
|5.2| their estimator shows a directional bias. We did not inves- 
tigate any further, whether this bias invalidates their findings at 
a rather low flux limit. Our analysis based on the three dimen- 
sional estimator applied to NVSS and using the mask defined by 
Singal|([20TT) gives (RA, dec) = (158° ± 21°, -2° ± 21°). 



Singal 2011 1 flux weighted number counts, SIF* corrects SIF 
for slope (this work), 3DS three-dimensional estimator, mask 



adopted from Singal (201 1 ) (this work), 2DCG two-dimensional 
estimator, CG mask (this work), GH ( Gibelyou & Huterer 2012 ); 
Radio dipole from WENSS: 2DCG two-dimensional estimator, 
CG mask (this work). The expectations for a purely kinetic radio 
dipole are given at the bottom of the table. 



Source 


Flux > 


N 


RA 


dec 


d 




(mJy) 




o 


o 


(io- 2 ) 


NVSS 


BW 


25 


197,998 


158 ±30 


-4 ±34 


1.1 ±0.3 


SIN 


25 


185,474 


158 ± 10 


-2+10 


1.8 ±0.4 


SIF 


25 


184, 237 


159 ± 10 


-7 ±9 


2.2 ±0.6 


SIF* 










1.6 ±0.5 


3DS 


25 


185,649 


158 ± 21 


-2 ±21 


1.6 ±0.6 


2DCG 


25 


195,245 


155 ± 14 




1.9 ±0.5 


GH 


15 


211,487 


117 ±20 


6± 14 


2.7 ±0.5 


3DS 


15 


298,289 


149 ± 19 


17 ±19 


1.4 ±0.5 


2DCG 


15 


313,724 


148 ± 15 




1.5 ±0.5 


WENSS 


2DCG 


25 


92,600 


117 ± 40 




2.9 ± 1.9 


expected 


NVSS 






168 


-7 


0.48 ± 0.04 


WENSS 






168 


-7 


0.42 ± 0.03 



Table 8. Comparison of results. Radio dipole from NVSS: BW 
Blake & WaIIl|2002] l, SIN ( |Singal|[20TTl l number counts, SIF 



For the amplitude of the radio dipole, the situation is more 
contrived. Here we focused on the study of linear estimators and 
showed that all linear estimators under investigation returned a 
biased estimate of the amplitude. The amplitude estimators of 
|Blake & Walrj ( p002] i and |Gibelyou & Huterer| ( |20T2] > are un- 
biased, but the latter one uses a biased direction estimate as an 
input and is thus of limited interest. Besides bias, we identified 
another effect that reduces the dipole amplitude found by the flux 
estimator used in Singal ( 201 1 1. We reduced the result of this es- 
timator by a factor of 1 .4, due to the fact that the appropriate ex- 
ponent of the differential number count is given by i = x + 1 (see 
section [6j. With this correction, the result of the flux weighted 
estimator (denoted by SIF* in table [8} is now in agreement with 
the result of |Blake & ' WaUlp002) . 

Our three dimensional estimate with the masking of Singal 



(2011 



gives rise to d — (1.6 ± 0.6) X 10" 2 . One should keep in 
mind that this is a biased result, thus one cannot naively compare 
it to the expected amplitude. In order to figure out if our result 
is consistent with the null hypothesis that the radio sky is sta- 
tistically isotropic, modified by the kinetic effects of our proper 
motion (measured via the CMB dipole), we performed 100,000 
Monte Carlo simulations. The corresponding histogram is shown 
in figure [10] We find that only 120 of those realizations contain 
a dipole higher than the measured one and thus we can exclude 
that the estimated radio dipole is just due to our proper motion 
and amplitude bias at 99.5% CL. This is actually very puzzling, 
as the direction of the radio dipole agrees with the direction of 
the CMB dipole within the measurement error. 

We can redo this analysis with the null hypothesis that the 
radio dipole was accurately measured by |Blake & Wall (2002). 
This time we find that 8836 out of 100,000 realisations are higher 
than our measured dipole. If we increase the implemented veloc- 
ity towards the upper one sigma bound of the dipole amplitude 
from Blake & Wall (2002), we observe every third simulation to 
be above our own measurement (30.0%). Therefore our result is 
in agreement with Blake & Wall (2002). 

Before we turn to the discussion of potential explanations, 
let us inspect the dipole amplitudes from the two dimensional 
estimator. The dipole amplitudes estimated with the two dimen- 



sional estimator are also in agreement with the results of |Singal| 
pOTTl ) and |Blake & Wall|p05^ . We find d sin § d = (1.9 + 0.5)x 
lF^for the NVSS analysis and dsin§ d = (2.9 + 1.9) x 10~ 2 
for WENSS, which translate into lower limits on d. Thus, the 
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results from the WENSS catalogue are in agreement with the 
radio dipole found in the NVSS catalogue. This is encouraging 
as they are prepared at different instruments and probe different 
frequencies. 

In figures [10] and [TTJ we present the corresponding results 
from 100,000 Monte Carlo simulations for the geometries of the 
NVSS and WENSS two dimensional estimators. In both cases 
we find that the amplitude bias is not enough to explain the dif- 
ference between the observed and the expected amplitude. In 
the case of NVSS the null hypothesis (isotropic sky plus proper 
motion dipole) is ruled out at 99.6% CL, while for the WENSS 
analysis it is inconsistent at 98.1% CL. 

The two dimensional estimation (2DCG), using the NVSS 
catalogue, is also compared to simulations, assuming the results 
from |Blake & Wall| | |2002) . We observe 1141 and 7024 out of 
100,000 simulations to give a dipole above our measurement for 
the assumptions of d tme = 1.1 x 10~ 2 and ii trae = 1.4 x 10~ 2 , 
respectively. Therefore, the result of our two dimensional esti- 
mator is not in contradiction to Blake & Wall ( 2002 ). 



0.06 



0.05 



0.04 



0.03 - 



3 0.02 



0.01 



simulation C 



0.005 0.01 0.015 0.02 0.025 0.03 0.035 0.04 
d 

Fig. 11. Histogram of dipole amplitudes from 100,000 simula- 
tions for the dimensional estimator, assuming the CMB expecta- 
tion, a slope of x = 0.8, 92,482 sources per simulation and the 
CG masking form for the WENSS catalogue. The black vertical 
line is the WENSS result of this work. 



9. Conclusion 

We conclude that the task to measure the cosmic radio dipole re- 
mains interesting and we expect that interesting information on 
cosmology can be extracted from this measurement. All mea- 
surements so far point towards a higher radio dipole amplitude 
than expected, when we assume that the cosmic radio dipole is 
just due to our peculiar motion with respect to the rest frame de- 
fined by the CMB. This is quite puzzling, as the orientation of 
the radio dipole agrees with the orientation of the CMB dipole 
within measurement errors. This is the case for the NVSS and the 
WENSS analysis, two radio point source catalogues that cover 
~ 3n and ~ n of the sky, respectively. They contain information 
at different frequencies (1.4 GHz and 325 MHz) and have been 
put together by different instruments and thus provide a strong 
constraint on systematic issues. 

Our detailed analysis of various linear dipole estimators 
( |Crawford||2009l |Singal|[20TT] |Gibelyou & Huterer|[20T2l i for 



three dimensional estimates (aj, 6j, d) and a new linear estima- 
tor for a two dimensional estimate {ad, d cos 5d) had to tackle 
several non-trivial issues. We investigated issues of directional 
bias, amplitude bias and masking. There is still room to opti- 
mize the masking of the galaxy. We did not look into quadratic 
estimators, as used by Blake & Wall (2002). Our studies did not 



incorporate the uncertainties in point source positions, as we ex- 
pect that they are subdominant (their magnitude is well below 
the effect of aberration, see e.g. |Rengelin k et al.| (fT997l). The 
measurement error of the flux is also expected to be subdominant 
as we included sources with flux above 15 mJy only. In the case 
of NVSS these are a factor 6 above the 5cr point source detection 
limit, in case of WENS S it is a factor of 1 .4. The number of point 
sources considered in our analysis is about 190,000 for NVSS 
and 92, 000 for WENSS. Putting all those facts together, we con- 
sider the NVSS analysis to be more reliable. Nevertheless, the 
results of the WENSS analysis are fully consistent with our re- 
sults from NVSS. 

Our result (3DS) for the NVSS catalogue is (RA, dec, d) = 
(154° + 21°, -2° ± 21°, (1.6 + 0.6) x 10" 2 ). Thus we conclude 
that the observed amplitude of the radio dipole exceeds the ex- 
pected amplitude by about a factor of three. We could imagine 
that this might be due to the structure that causes our proper mo- 
tion, which in a simple model of our Hubble patch would cer- 
tainly be aligned with the direction of proper motion. However, 
all attempts to identify a "great attractor" by means of other ob- 
servations (optical, infra-red, X-ray) have failed so far. Of all 
those probes, the radio surveys are definitely the deepest probe 
of the Universe, as the mean redshift of NVSS sources is esti- 



mated to be 1.2 by |de Zotti et al.| ( [20T0[ i and 1.5 by |Ho et aL 
(2008). In order to explain the observed excess radio dipole by 
contributions from local structure, we would need a density con- 
trast of order 0.05 at scales extending to about z ~ 0.3, which 
does not seem plausible. Without a detailed study of the redshift 
distribution of the radio sources it is impossible to judge whether 
this finding is actually in agreement with the current standard 
model of cosmology. 

Another reason for the large amplitude of the radio dipole 
could be that the linear estimators considered in this work do 
not assume the deviation from isotropy to be a pure dipole. Thus 
higher multipole moments also contribute to the measured am- 
plitude. 

It has been found from the analysis of the CMB that 
quadrupole, octopole and a few more low ^-multipoles seem not 
to be orientated randomly on the sky, but show some unexpected 
alignments ( |Schwarz et al.|2004 Bennett et al.|201 1" Copi et al. 



2010) among themselves and with the CMB dipole direction. 



Thus it might not be surprising that also the dominant anisotropy 
direction of the radio sky lines up. 

It is evident that it would be crucial to repeat this investiga- 
tion with new and even deeper radio catalogues, which provide 
more sources. In the near future there will be three large sky 
surveys available ( jRaccanelli et al.||2012| ). A multi-wavelength 
study will be possible based on the Multifrequency Snapshot 
Sky Survery (MSSS) of the International LOFAR Telescope and 
with the LOFAR Tier-1 survey. The Australian Square Kilometre 
Array (SKA) Pathfinder (ASKAP) will produce the Evolutionary 
Map of the Universe (EMU) and the Westerbrok Synthesis Radio 
Telescope (WSRT) equipped with the Apertif system will com- 
pile the Westerbrok Observations of the Deep Apertif Northern 
sky survey (WODAN) catalogue. 

The multi-wavelength surveys will also allow us to directly 
measure the spectral index a, which has to be known to connect 
the measured amplitude to the kinetic dipole. A steepening of 
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the spectral index for the lowest flux sources would increase the 
expected kinetic amplitude. HI surveys will have the advantage 
that they will also provide redshift information on top of posi- 
tions and fluxes and we will be able to study the evolution of the 
radio dipole as a function of redshift out to redshifts of a few. 
Beyond that SKA will increase the number of sources in such a 
survey by orders of magnitude. All these surveys will reduce the 
random dipole contribution, improve on systematics, and allow 
us to settle the question: Is the radio dipole in agreement with 
the CMB dipole? 
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Appendix A: Monte Carlo Simulations 

We used the pseudorandom number generator Mersenne Twister 
for all Monte Carlo simulations. One simulated source consists 
of two position coordinates and one flux value. The coordi- 
nates will be drawn from an uniform distribution, leading to 
an isotropic sky for a catalogue with many sources. In order to 
obtain a desired number count power law like Q with a cer- 
tain slope x, we calculate the flux S using a random number A 
(choosen between and 1) by 



S L (l-Ay 



(A.l) 



where Si is a flux value 20% below the simulated flux limit. 
The simulation also creates sources, which only due to Doppler 
shifting are counted in the final catalogue, because this value 
of Si lies below the simulated flux limit. Lowering Si further 
increases computational time and is not necessary, as long as the 
simulated velocity v is below 0.1c. 

The two physical effects (Doppler shift and spherical aberra- 
tion) will be implemented separately. In cooperating the Doppler 
effect is straightforward, since it only affects the flux values of 
each source depending on the angle 9 between the source and 
the velocity direction, i.e. 



S obs(Vobs) = S r 



1 + v - cos(6 l ) 



7 



v\2 



(A.2) 



In the simulation a is fixed to 0.75 for all sources. The velocity 
direction and amplitude can by chosen be the user. 

To model the relativistic effect of stellar aberration one has to 
change the position of each radio source. The aberration formula 
is 



tan(f/) 



sin(0)^ 



1-0 



r ~v2 



+ cos(fT) 



(A.3) 



where 9' is the new angle between the velocity direction and 
the radio source. Forth, the position of a radio source is trans- 
lated into Cartesian coordinates by assuming that it lies on a unit 
sphere. Then a straight line from this point (P) to the velocity 



direction on the sphere (V) is constructed depending on a pa- 
rameter t 



r(t) = P(l-i) + Vt. 



(A.4) 



On this line we choose a f in such a way that r(t') points towards 
the new position. The value of t' can be determined by 



r(t') ■ V = r'(t') cos(6»') 
, _ r'(t') cos(f/) - cos(f/) 
~* ' ~ 1 - cos(fT) 

with r'(f') = ^r' 2 (f). This equation is solved by 
sin(6»-f/) , sin(0 + 0') 



sin(6> - 0') + sin(6>') 



sin(0 + 0') -sin(6»') 



(A.5) 
(A.6) 



(A.7) 



We know that for 9-9' the result of t' must always be 0. 
Therefore the correct solution is t' - t\. Now one has to trans- 
form r(t') back into spherical coordinates in order to find the new 
position of the radio source. The new declination is then (the 
index v stands for the velocity direction): 

=arccos^[(l-f')cos(i?) + f'cos(# v )]j (A.8) 



and the new right ascension ip': 

'(1 - t') sin(??) sin(<p) + f sin(# v ) sin(a v ) 



ip = arcstn 



r' sin(??') 



(A.9) 



This way one obtains a simulated sky, including the effect of the 
observers movement. Now one can feed the different estimators 
with those sky simulations and obtain the resulting dipole vec- 
tors. 
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